function [outdata] = compare_plot(varargin)
%
% Function to make plots that compare data across entries or cells
%
persistent rtitle result state_arg

global CONTROL
global TRES EPSC_STATE

outdata=[];
msize = 3;
fsize = 7;
if(nargin > 0)
    arg = varargin{1};
else
    arg= '';
end;

if(nargin > 1)
    arg2 = varargin{2};
else
    arg2 = 0;
end;

if (isempty(state_arg))
    state_arg = 0;
end;
EPSC_STATE=state_arg;

h2=findobj('Tag', 'CTL'); % first get the selection number
sel = get(h2, 'Value');
if(length(sel) < 1) % must have something to combine in selection...
    msgbox('Select Data Sets First', 'Error!', 'replace');
    return;
end;
% make figure for plots
hf = findobj('Tag', 'Comparison');
if(isempty(hf))
    hf = figure('Tag', 'Comparison', 'Name', 'Comparison Plots', 'NumberTitle', 'off');
end;
figure(hf);
tf=findobj('tag', 'Compare_export');
if(~strcmp(arg, 'export') & isempty(tf)) % add an export data command...
    %clf;
    f = uimenu('Label', 'E&xport Data', 'Position', 6, 'tag', 'Compare_export');
    uimenu(f, 'Label', '&Export', 'Callback', 'compare_plot export;');
end;
leg = '';
% fprintf(1, 'arg = %s\n', arg);
switch (arg)

    case 'clear' % clear the comparison list results...
        result = [];
        rtitle = [];
        clf;

    case 'plot_result'
        TRES = result;
        s = size(result);
        tb = result(1,:);
    case 'return_result'
        outdata = result;
    

    case('iv') % multiply selected IVs in one graph
        l = 1;
        leg={};
        result ={};
        rtitle = 'IV';
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'iv') & ~isempty(CONTROL(j).iv))
                plot(CONTROL(j).iv.mean(1,:), CONTROL(j).iv.mean(2,:), 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                leg{l} = sprintf('%s %s',CONTROL(j).filename, CONTROL(j).reclist);
                l = l + 1;
                hold on;
                if(j == sel(1))
                    result = CONTROL(j).iv.mean(1,:);
                end;

                result = {result; CONTROL(j).iv.mean(2,:)};
            end;
        end;
        hold off;
        legend(leg,2);

    case('fi') % multiply selected spike rate vs current (FI) plot
        l=1;
        leg={};
        result ={};
        for j=sel % copy the information to the main record
            ks = l; %j - sel(1) + 1;
            s = symbols(-ks);
            if(check_field(CONTROL(j), 'iv') && ~isempty(CONTROL(j).iv))
                plot(CONTROL(j).iv.mean(1,:), CONTROL(j).iv.nr_spk, 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                leg{l} = sprintf('%s : %s',CONTROL(j).filename, CONTROL(j).reclist);
                l = l + 1;
                hold on;
                if(j == sel(1))
                    result = CONTROL(j).iv.mean(1,:);
                end;

                result = [result; CONTROL(j).iv.nr_spk];
            end;
        end;
        hold off;
        legend(leg,2);

    case 'cciv'
        mksize = 3;
        l = 1;
        subplot(2,2,1);
        % first compare the IVs
        result = [];
        rtitle = [];
        if(~isfield(CONTROL, 'iv')) % verify we have a field
            return;
        end;
        civ = [CONTROL(sel).iv]; % just get the stuff we need...
        % reduce selection to just those that really have an analysis in them.
        sels = [];
        sel2 = [];
        k=1;
        for j = sel
            if(isempty(CONTROL(j).iv))
                continue;
            else
                sel2 = [sel2 k];
                k = k + 1;
                sels = [sels j];
            end;
        end;
        exp = strvcat(CONTROL(getmainselection).E_C); % experimental conditions
        g1 =strvcat(CONTROL(getmainselection).G1);
        g2 = strvcat(CONTROL(getmainselection).G2);
        uexp = unique(cellstr(strvcat(exp, g1, g2)));

        for j=1:length(sel2) % copy the information to the main record
            ks = l; % l*2 - 1; % j - sel(1) + 1;
            if(length(uexp) < 2)
                s = symbols(-ks);
            else
                s = symbols(-j);
            end;
            plot(civ(j).mean(1,:), civ(j).mean(2,:), 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', mksize, 'LineStyle', '-');
            leg{l} = sprintf('%s %s',CONTROL(sels(j)).filename, CONTROL(sels(j)).reclist);
            l = l + 1;
            hold on;
        end;
        hold off;
        %  legend(leg,2);
        xlabel('Current (nA)', 'Fontsize', 8);
        ylabel('V (mV)', 'Fontsize', 8);
        title('IV');
        subplot(2,2,2);
        % next compare the FI plots
        l = 1;
        for j=1:length(sel2) % copy the information to the main record
            ks = l; %l*2 - 1; % j - sel(1) + 1;
            s = symbols(-ks);
            plot(civ(j).mean(1,:), civ(j).nr_spk, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-', 'MarkerSize', mksize);
            %            leg(l,:) = sprintf('%s %s',CONTROL(sels(j)).filename, CONTROL(sels(j)).reclist);
            l = l + 1;
            hold on;
        end;
        hold off;
        legend(leg,2);
        xlabel('Current (nA)', 'Fontsize', 8);
        ylabel('Spike Count (N)', 'Fontsize', 8);
        title('FI');

        % fsl
        % first we have to get the data without the spikes.
        subplot(2,2,3);
        l = 1;
        for j=1:length(sel2) % copy the information to the main record
            dj = find(civ(j).nr_spk > 0);
            ks = l; % l*2 - 1; % j - sel(1) + 1;
            s = symbols(-ks);
            plot(civ(j).mean(1, dj), civ(j).fsl, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', mksize, 'LineStyle', '-');
            %           leg(l,:) = sprintf('%s %s',CONTROL(sels(j)).filename, CONTROL(sels(j)).reclist);
            l = l + 1;
            hold on;
        end;
        hold off;
        %   legend(leg,2);
        xlabel('Current (nA)', 'Fontsize', 8);
        ylabel('FSL(ms)', 'Fontsize', 8);
        title('FSL');

        % fisi
        subplot(2,2,4);
        l = 1;
        for j=1:length(sel2) % copy the information to the main record
            dj = find(civ(j).nr_spk > 1);
            ks = l; %l*2 - 1; % j - sel(1) + 1;
            s = symbols(-ks);
            plot(civ(j).mean(1, dj), civ(j).fisi, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', mksize, 'LineStyle', '-');
            %          leg(l,:) = sprintf('%s %s',CONTROL(sels(j)).filename, CONTROL(sels(j)).reclist);
            l = l + 1;
            hold on;
        end;
        hold off;
        %   legend(leg,2);
        xlabel('Current (nA)', 'Fontsize', 8);
        ylabel('FISI(ms)', 'Fontsize', 8);
        title('FISI');

    case('averageiv')
        % plot the average IV and the standard deviation
        % Frst we interpolate all the voltages onto a commmon current scale
        % then we do the average.
        rtitle = 'AverageIV';
    exp=[]; g1=[]; g2 = [];
        civ = [CONTROL(getmainselection).iv];
%        exp = strvcat(CONTROL(getmainselection).E_C); % experimental conditions
%        g1 =strvcat(CONTROL(getmainselection).G1);
        g2 = strvcat(CONTROL(getmainselection).G2);
        uexp = unique(cellstr(strvcat(exp, g1, g2)));
        for k = 1:length(uexp)
            expsel = strmatch(uexp(k), g2);
            if(isempty(expsel))
                continue;
            end;
            mxi = 0;
            mni = 0;
            np = zeros(length(expsel), 1);
            for j = 1:length(expsel)
                i = expsel(j);
                mi = max(civ(i).mean(1,:));
                if(mi > mxi)
                    mxi = mi;
                end;
                mi = min(civ(i).mean(1,:)); % get max and mi of currents
                if(mi < mni)
                    mni = mi;
                end;
                np(j) = length(civ(i).mean(1,:));
            end;
            mnp = ceil(mean(np)); % get the mean number of points
            increment = 20; % incmement in pA
            mxi = increment*ceil(mxi/increment);
            mni = increment*floor(mni/increment);
            Iref = mni:increment:mxi;
            V = zeros(length(expsel), length(Iref));
            for j = 1:length(expsel)
                i = expsel(j);
                V(j,:) = interp1(civ(i).mean(1,:), civ(i).mean(2,:), Iref);
            end;
            Vavg = zeros(size(V,2), 1);
            var = zeros(size(V,2), 1);
            navg = var;
            % adjust relative to resting membrane potentail
            [x,  jRm] = min(abs(Iref));
            Vrm = V(:,jRm);
            for i = 1: size(V, 1)
                V(i,:) = V(i,:) - Vrm(i);
            end;
            for i = 1: size(V,2)
                [Vavg(i), var(i), navg(i)]  = mean_var(V(:,i));
            end;


            include = find(navg >= 3); % fewer than 3 observations at any point - exclude it
            Vavg = Vavg(include);
            var = var(include);
            Iref = Iref(include);
            s = symbols(k);
            Vstd = sqrt(var);
            h2 = plot(Iref, Vavg);
            set(h2, 'color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
            hold on;
            h = my_errorbar(Iref, Vavg, -Vstd, Vstd);
            set(h, 'color', s.c);
            hold on;
            if(arg2 > 0)
                for i = 1:length(Iref)
                    fprintf(1, '%f\t%f\t%f\n', Iref(i), Vavg(i), Vstd(i));
                end;
            end;
        end;


    case('averagespikecount')
        % plot the average IV and the standard deviation
        % Frst we interpolate all the voltages onto a commmon current scale
        % then we do the average.
        % note: for this we eliminate the non-monotonic parts of the curve
        % where depolarization drives the cell into Na channel block

        rtitle = 'AverageSpikeCount';
 %       clf;
        hc =  multiplot('Comparison', 3, 1, 1);
        civ = {CONTROL(getmainselection).iv};
    exp=[]; g1=[]; g2 = [];
%        exp = strvcat(CONTROL(getmainselection).E_C); % experimental conditions
%        g1 =strvcat(CONTROL(getmainselection).G1);
        g2 = strvcat(CONTROL(getmainselection).G2);
        uexp = unique(cellstr(strvcat(exp, g1, g2)));
        for k = 1:length(uexp)
            expsel = strmatch(uexp(k), g2);
            if(isempty(expsel))
                continue;
            end;
            mxi = 0;
            mni = 0;
            np = zeros(length(expsel), 1);
            for j = 1:length(expsel)
                i = expsel(j);
                   mi = max(cell2mat({civ{i}.mean(1,:)}));
                if(mi > mxi)
                    mxi = mi;
                end;
                    mi = min(cell2mat({civ{i}.mean(1,:)}));
                if(mi < mni)
                    mni = mi;
                end;
                np(j) = size(civ{i}.mean, 2);
            end;
            mnp = ceil(mean(np)); % get the mean number of points
            increment = 20; % incmement in pA
            mxi = increment*ceil(mxi/increment);
            mni = increment*floor(mni/increment);
            Iref = mni:increment:mxi;
            V = zeros(length(expsel), length(Iref));
            for j = 1:length(expsel)
                i = expsel(j);
                V(j,:) = interp1(cell2mat({civ{i}.mean(1,:)}), cell2mat({civ{i}.nr_spk}), Iref);
            end;
            maxi = size(V, 2);
            if(k == 1)
                maxnspk = zeros(size(V, 1), 1);
            end;
            for i = 1: size(V, 1)
                [maxn, pnspk] = max(V(i, :));
                if(pnspk < maxi)
                    V(i, pnspk+1:end) = NaN;
                end;
                if(k == 1)
                    maxnspk(i) = maxn;
                end;
                V(i,:) = V(i,:)/maxnspk(i); % normalize by max spike count

            end;
            V
            Vavg = zeros(size(V,2), 1);
            var = zeros(size(V,2), 1);
            navg = var;
            for i = 1: size(V,2)
                [Vavg(i), var(i), navg(i)]  = mean_var(V(:,i));
            end;
            include = find(navg >= 5) % fewer than 3 observations at any point - exclude it
            Vavg = Vavg(include);
            var = var(include);
            Iref2 = Iref(include);
            s2 = symbols(k);
            Vstd = sqrt(var);
            subplot(hc(1));
            hold on;
            h2 = plot(Iref2, Vavg);
            Iref2
            Vavg
            Vstd
            h = my_errorbar(Iref2, Vavg, -Vstd, Vstd);
            set(h2, 'color', s2.c, ...
                'Marker', s2.s, 'MarkerFaceColor', s2.f, 'LineStyle', '-');
            set(h, 'color', s2.c, 'LineStyle', '-');
            u = get(gca, 'Xlim');
            set(gca, 'Xlim', [0 u(2)]);
            subplot(hc(3));
            for j = 1: length(expsel)
                s = symbols(j);
                i = expsel(j);
                plot(Iref, V(j,:), 'color', s2.c, 'MarkerSize', 8, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                hold on;
            end;
            u = get(gca, 'Xlim');
            set(gca, 'Xlim', [0 u(2)]);
            VS{k} = V;
            IS{k} = Iref;
        if(arg2 > 0)
            for i = 1:length(Iref2)
                fprintf(1, '%f\t%f\t%f\n', Iref2(i), Vavg(i), Vstd(i));
            end;
        end;
        end;
        % lastly, a simple difference plot
        subplot(hc(2));
        vc= VS{1};
        vs = VS{2};
        vd = vs - vc;
        for j = 1: length(expsel)
            s = symbols(j);
            i = expsel(j);
            plot(IS{1}, vd(j,:), 'color', s2.c, 'MarkerSize', 8, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
            hold on;
        end;
        u = get(gca, 'Xlim');
        set(gca, 'Xlim', [0 u(2)]);


    case('CCIVSummary')
  %      clf;
        normalized = 0; % doesn't work yet...
        sf=getmainselection;
        rtitle = 'CCIV_Summary';
%        hc =  multiplot(rtitle, 1, 3, 1);
%        clf;
        hc =  multiplot(rtitle, 1, 3, 1);
        sf = getmainselection;
        civ = [CONTROL(sf).iv];
        exp = strvcat(CONTROL(sf).E_C); % treatments conditions
        uexp = unique(cellstr(strvcat(exp)));
        cond = strvcat(CONTROL(sf).G2);
        ucond = unique(cellstr(strvcat(cond))); % animal condition
        Rin_m = zeros(length(uexp), length(ucond));
        Rin_s = Rin_m; Rin_N = Rin_m;
        Rmp_m = Rin_m; Rmp_s = Rin_m; Rmp_N = Rin_m;
        Tau_m = Rin_m; Tau_s = Rin_m; Tau_N = Rin_m;
        for j = 1:length(ucond)
            condsel = strmatch(ucond(j), cond);
            for k = 1:length(uexp)
                expsel = strmatch(uexp(k), exp);
                if(isempty(expsel))
                    continue;
                end;
                sfsel = sf(intersect(expsel, condsel));
                if(normalized) % each cell serves as it's own control here
                    fns = {CONTROL(sfsel).filename};
                    fnl = unique(fns);
                    Rin_l = zeros(length(fnl), 1);
                    Rmp_l = Rin_l;
                    Tau_l = Rin_l;
                    for i = 1:length(fnl)
                        bma = strcmpi(fnl(i), {CONTROL.filename}); % find each position
                        ma = find(bma == 1);
                        sfsel2 = intersect(ma, sfsel);
                        Rin_l(i) = [CONTROL(sfsel2).Rin] / CONTROL(sfsel2(1)).Rin;
                        Rmp_l(i) = [CONTROL(sfsel2).Rmp] - CONTROL(sfsel2(1)).Rmp;
                        Tau_l(i) = [CONTROL(sfsel2).iv.tau1]/ CONTROL(sfsel2(1)).iv.tau1;
                    end;

                    [Rin_m(k, j), Rin_s(k, j), Rin_N(k, j)] = mean_var(Rin_l);
                    [Rmp_m(k, j), Rmp_s(k, j), Rmp_N(k, j)] = mean_var(Rmp_l);
                    [Tau_m(k, j), Tau_s(k, j), Tau_N(k, j)] = mean_var(Tau_l);

                else
                    [Rin_m(k, j), Rin_s(k, j), Rin_N(k, j)] = mean_var([CONTROL(sfsel).Rin]);
                    [Rmp_m(k, j), Rmp_s(k, j), Rmp_N(k, j)] = mean_var([CONTROL(sfsel).Rmp]);
                    [Tau_m(k, j), Tau_s(k, j), Tau_N(k, j)] = mean_var([civ(expsel).tau1]);
                end;
            end;
        end;
        XE = repmat([1:length(uexp)]', 1, length(ucond));
        X = XE(:,1);
        BW = 0.8;


        Rin_s = sqrt(Rin_s);
        Rmp_s = sqrt(Rmp_s);
        Tau_s = sqrt(Tau_s);
        subplot(hc(1));
        hb = bar(X, Rin_m, BW);
        hbc = get(hb, 'Children');
        if(length(hbc) == 1)
            hbc = {hbc};
        end;
        for i = 1: size(XE, 1)
            for j = 1: size(hbc)
                xd = get(hbc{j}, 'XData');
                XE(i,j) = deal(mean(xd(:,i))); % use middle of plotted bar
            end;
            hold on;
        end;
        he = my_errorbar(XE, Rin_m, -Rin_s, Rin_s, '.');
        ht=title('Rin (Mohm)');
        py = get(ht, 'Position');
        set(ht, 'VerticalAlignment', 'top', 'Position', [py(1) 0.95*py(2) py(3)]);
        set(gca, 'xticklabel', uexp);
        subplot(hc(2));
        bar(X, Rmp_m, BW);
        hold on;
        he = my_errorbar(XE, Rmp_m, -Rmp_s, Rmp_s, '.');
        ht=title('RMP (mV)');
        py = get(ht, 'Position');
        set(ht, 'VerticalAlignment', 'top', 'Position', [py(1) 0.95*py(2) py(3)]);

        set(gca, 'xticklabel', uexp);
        subplot(hc(3));
        bar(X, Tau_m, BW);
        hold on;
        he = my_errorbar(XE, Tau_m, -Tau_s, Tau_s, '.');
        ht=title('Taum (ms)');
        py = get(ht, 'Position');
        set(ht, 'VerticalAlignment', 'top', 'Position', [py(1) 0.95*py(2) py(3)]);
        set(gca, 'xticklabel', uexp);
        drawnow;


    case('rawact') % multiply selected IVs in one graph
        l = 1;
        rtitle = 'RawAct';
        for j=sel % copy the information to the main record
            ks = l; % (j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'act') && ~isempty(CONTROL(j).act))
                plot(CONTROL(j).act.ss_V, CONTROL(j).act.ss_I, 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                hold on;
                leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
                %      s = symbols(ks+1);
                %      plot(CONTROL(j).act.pk_V, CONTROL(j).act.pk_I, 'Color', s.c, ...
                %          'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                %      leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                %      l = l + 1;
                if(j == sel(1))
                    result = CONTROL(j).act.ss_V;
                end;

                result = [result; CONTROL(j).act.ss_I];
            end;
        end;
        legend(leg);
        hold off;

    case ('act') % normalized activation to +20 mV
        l = 1;
        rtitle = 'Normalized Activation';
        normv = 40; % test normalization voltage
        testv = -50; % voltage to measure test fraction of activation at
        vpct = [];
        vfn = [];
        Vr = -85;
        for j=sel % copy the information to the main record
            ks = l; % (j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'act') & ~isempty(CONTROL(j).act))
                IO = CONTROL(j).act.pk_V;
                G = CONTROL(j).act.pk_I;  % Don't transform to conductance
                [x, ind] = min(abs(IO-normv)); % find point closest to normalizing voltage
                if(ind < length(CONTROL(j).act.pk_I))
                    indxs=[ind-1:ind+1];
                else
                    indxs=[ind-1:ind];
                end;
                Gmax = mean(G(indxs));
                iok=find(abs(IO-Vr)>=10);
                %  IO = IO(iok);
                %   G = G(iok);
                hold on;
                plot(IO, G/Gmax, 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', 'none');
                leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                [y, yind] = min(abs(IO - testv));
                % vpct(l) = {CONTROL(j).filename, G(yind)/Gmax};
                vfn(l) = j;
                vpct(l) =  G(yind)/Gmax;
                l = l + 1;
                % fit the data in the plot to a double boltzmann...
                a0 = 0; a1 = 0.5; k1 = 2; v1 = -84; x1  = -60; k2=2; x2=10;
                [ssiv_par, ssiv_chisq, ssiv_niter, ssiv_fit] = mrqfit('boltzmann', [a0 a1 k1 x1 v1 ], IO, G/Gmax, [], ...
                    [0 1 1 1 0], [], [], 500, []);
                plot(IO, ssiv_fit, 'linestyle', '-','Color', s.c, ...
                    'Marker', 'none');
                %ssiv_par
            end;
        end;
        set(gca, 'YLim', [0 1.2]);
        legend(leg,2);
        hold off;
        for i = 1:length(vpct)
            fprintf(1, '%s \t %6.3f\n', CONTROL(vfn(i)).filename, vpct(i))
        end;
        % result = [vfn'; vpct'];
        % rtitle = ['filename'; sprintf('Activation at %d mV', testv )];

        % list the 20% of max voltage.


    case('inact') % multiply selected IVs in one graph
        l = 1;
        rtitle = 'Inactivation';
        for j=sel % copy the information to the main record
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'inact') & ~isempty(CONTROL(j).inact))
                %plot(CONTROL(j).inact.ss_V, CONTROL(j).inact.ss_I/max(CONTROL(j).inact.pk_I), 'Color', s.c, ...
                %  'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                hold on;
                leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
                s = symbols(ks+1);
                if(min(CONTROL(j).inact.pk_I) > 0)
                    plot(CONTROL(j).inact.pk_V, CONTROL(j).inact.pk_I/max(CONTROL(j).inact.pk_I), 'Color', s.c, ...
                        'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                    leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                    l = l + 1;
                end;
            end;
        end;
        hold off;
        legend(leg,1);

    case('spike') % multiple selected averaged spikes in one graph
        l = 1;
        result = [];
        rtitle = 'spike';
        for j=sel % copy the information to the main record
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'spike') && ~isempty(CONTROL(j).spike)) % spikes are plotted with no symbols
                plot(CONTROL(j).spike.tavg, CONTROL(j).spike.vavg, 'Color', s.c);
                hold on;
                leg{l} = sprintf('%12s %s',CONTROL(j).reclist, CONTROL(j).Commentary);
                l = l + 1;
                s = symbols(ks+1);
                if(j == sel(1))
                    %    result = CONTROL(j).spike.tavg;
                end;
                %result = [result; CONTROL(j).spike.vavg];
            end;
        end;
        hold off;
        legend(leg,1);

    case('ispike') % multiple selected averaged current preceeding a spike in one graph
        l = 1;
        result = [];
        rtitle = 'ispike';
        for j=sel % copy the information to the main record
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'spike') & ~isempty(CONTROL(j).spike)) % spikes are plotted with no symbols
                plot(CONTROL(j).spike.tavg, CONTROL(j).spike.iavg, 'Color', s.c);
                hold on;
                leg{l} = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
                s = symbols(ks+1);
                if(j == sel(1))
                    result = CONTROL(j).spike.tavg;
                end;
                result = [result; CONTROL(j).spike.iavg];
            end;
        end;
        hold off;
        legend(leg,1);

    case {'hyp'}
        subplot('Position', [0.1,0.55,0.6,0.40]);
        l = 1;
        lg = cell(length(sel),1);
        result = {};
        rtitle = 'hyp';
        for j=sel % copy the information to the main record
            k = l; % find(j == sel);
            s = symbols(-k);
            if(check_field(CONTROL(j), 'iv') & ~isempty(CONTROL(j).iv) & check_field(CONTROL(j).iv, 'mv_fsl'))
                hold on
                plot(CONTROL(j).iv.mv_fsl, CONTROL(j).iv.fsl, 'color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-', ...
                    'MarkerSize', msize);
                if(check_field(CONTROL(j).iv,  'hypb_Vh') & ~isnan(CONTROL(j).iv.hypb_Vh))
                    %   plot(CONTROL(j).iv.hypb_vfit, CONTROL(j).iv.hypb_lfit, '-r');
                end;
                set(gca, 'FontSize', fsize);
                %  grid;
                ylabel(sprintf('FSL (%s)', CONTROL(j).T_Unit));
                xlabel(CONTROL(j).V_Unit);
                set(gca, 'FontSize', fsize);
            else
                plot(-110, 20, 'color', s.c, 'Marker', s.s, 'MarkerFaceColor', s.f, 'Markersize', msize);
            end;
            bn = sprintf('%s-%d',CONTROL(j).filename, CONTROL(j).recbeg);
            rtitle = [rtitle sprintf('%sV, %sL, ', bn, bn)];
            result{k} = [CONTROL(j).iv.mv_fsl; CONTROL(j).iv.fsl];
            cm = CONTROL(j).Commentary;
            lg{k}=sprintf('%2d  %s ',k, char([CONTROL(j).Commentary]'));
            l = l + 1;
        end;
        hold off;
        legend(char(lg),0);
        set(gca, 'Xlim', [-100 -40]);
        subplot('Position', [0.1,0.05,0.6,0.40]);
        l = 1;
        for j=sel % copy the information to the main record
            s = symbols(-l);
            if(check_field(CONTROL(j), 'iv') & ~isempty(CONTROL(j).iv) & check_field(CONTROL(j).iv, 'mv_fisi'))
                hold on
                plot(CONTROL(j).iv.mv_fisi, CONTROL(j).iv.fisi, 'color', s.c, ...
                    'Marker', s.s, 'MarkerfaceColor', s.f, 'LineStyle', '-', ...
                    'MarkerSize', msize);
                set(gca, 'FontSize', fsize);
                %  grid;
                ylabel(sprintf('FISI (%s)', CONTROL(j).T_Unit));
                xlabel(CONTROL(j).V_Unit);
                set(gca, 'FontSize', fsize);
            else
                plot(-110, 20, 'color', s.c, 'Marker', s.s, 'MarkerFaceColor', s.f, 'Markersize', msize);
            end;
            lg{k}=sprintf('%2d: %c',k, char([CONTROL(j).Commentary]'));
            l = l + 1;
        end;
        hold off;
        set(gca, 'Xlim', [-100 -40]);
        % information...
        subplot('Position', [0.73, 0.05, 0.20, 0.97]);
        axis([0,1,0,1])
        axis('off')
        lys=0.020;
        yp=1;
        l = 1;
        txts = 6;
        for  j = sel
            k = find(j == sel);
            text(0,yp-l*lys,sprintf('%3d:  %s',k, CONTROL(j).filename),'FontSize', txts);
            %get the time from the notefile
            m = [];
            if(~isempty(CONTROL(j).NOTEFILE))
                m = find(number_arg(CONTROL(j).recbeg) == [CONTROL(j).NOTEFILE.frec]);
            end;
            if(isempty(m))
                tm = '??:??:??';
            else
                tm  = CONTROL(j).NOTEFILE(m).etime;
            end;
            text(0,yp-(l+1)*lys,sprintf('      [%4d-%4d] at %8s', ...
                number_arg(CONTROL(j).recbeg),number_arg(CONTROL(j).recend),tm), 'FontSize', txts);
            l = l + 2;

        end


    case('spect') % multiple spectra in one graph
        l = 1;
        pc = 0.65; % patch color (uniform grey)
        for j=sel % First plot the backgrounds (SEMs) if they exist...
            if(check_field(CONTROL(j), 'spike') & ~isempty(CONTROL(j).spike)) % spikes are plotted with no symbols
                if(check_field(CONTROL(j).spike, 'spect_sem'))
                    [f, p] =make_patch(CONTROL(j).spike.fspect, CONTROL(j).spike.spect, CONTROL(j).spike.spect_sem);
                    patch(f, p, [pc pc pc], 'EdgeColor', 'none');
                    hold on;
                end;
            end;
        end;
        for j=sel % Next Plot the lines on top of the variances...
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'spike') & ~isempty(CONTROL(j).spike)) % spikes are plotted with no symbols
                plot(CONTROL(j).spike.fspect, (CONTROL(j).spike.spect), 'Color', s.c);
                hold on;
                leg{l} = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
            end;
        end;
        hold off;
        legend(leg,1);

    case 'vm-specamp' % multiple vm - spectrum amplitude plots on one graph.
        l = 1;
        for j=sel % Next Plot the lines on top of the variances...
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'spike') & ~isempty(CONTROL(j).spike)) % spikes are plotted with no symbols
                subplot(2,1,2);
                plot(CONTROL(j).spike.pre_vm, CONTROL(j).spike.pre_specamp,'Color', s.c, 'marker', s.s, 'LineStyle', 'none');
                hold on;
                set(gca, 'YLim', [0 1.0]);
                set(gca, 'XLim', [-60 -35]);

                subplot(2,1,1);
                plot(CONTROL(j).spike.pre_vm, CONTROL(j).spike.pre_slope, 'Color', s.c, 'marker', s.s, 'LineStyle', 'none');
                set(gca, 'YLim', [0 0.25]);
                set(gca, 'XLim', [-60 -35]);
                hold on;
                leg{l} = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
            end;
        end;
        hold off;
        legend(leg,1);

    case 'SAC' % multiple SAC - shuffled autocorrelation analysis in one plot for comparison
        l = 1;
        result = [];
        rtitle = 'SAC';

        for j=sel % Next Plot the lines on top of the variances...
            ks = l; % 2*(j - sel(1)) - 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'SAC') & ~isempty(CONTROL(j).SAC)) % spikes are plotted with no symbols
                subplot(2,1,2);
                hp = semilogx(CONTROL(j).SAC.hx1, CONTROL(j).SAC.hy1);
                set(hp,'Color', s.c, 'marker', s.s, 'LineStyle', 'none');
                hold on;
                yg1 = gaussfunc(CONTROL(j).SAC.hx1, CONTROL(j).SAC.Gfit1);
                semilogx(CONTROL(j).SAC.hx1, yg1);

                %set(gca, 'XLim', [-60 -35]);

                subplot(2,1,1);
                hp = semilogx(CONTROL(j).SAC.hx2, CONTROL(j).SAC.hy2);
                set(hp, 'Color', s.c, 'marker', s.s, 'LineStyle', 'none');
                %set(gca, 'YLim', [0 0.25]);
                %set(gca, 'XLim', [-60 -35]);
                hold on;
                yg2 = gaussfunc(CONTROL(j).SAC.hx2, CONTROL(j).SAC.Gfit2);
                semilogx(CONTROL(j).SAC.hx2, yg2);
                leg{l} = sprintf('%12s',CONTROL(j).reclist);
                if(j == sel(1))
                    result = CONTROL(j).SAC.hx1;
                end;
                result = [result; CONTROL(j).SAC.hy1'; yg1; CONTROL(j).SAC.hy2'; yg2];

                l = l + 1;
            end;
        end;
        hold off;
        legend(leg,1);


    case 'tc_int'
        % multiple tc_ints on one graph....... in time... UPDATED 6/15/2006
        % P. Manis (Major changes).
        if(~isfield(CONTROL, 'tc_int'))
            return;
        end;
        % build plot window
        hp1 = subplot('Position', [0.1, 0.80, 0.8, 0.15]); % fsl
        hp2 = subplot('Position', [0.1, 0.63, 0.8, 0.15]); % fisi
        hp3 = subplot('Position', [0.1, 0.46, 0.8, 0.15]); % nspikes
        hp4 = subplot('Position', [0.1, 0.34, 0.8, 0.10]); % vm
        hp5 = subplot('Position', [0.1, 0.22, 0.8, 0.10]); % ihold
        hp6 = subplot('Position', [0.1, 0.10, 0.8, 0.10]); % rin

        l = 1;

        fns = unique({CONTROL(sel).filename}); % find out how many unique files are represented
        for k = 1:length(fns)
            sel2 = strmatch(fns{k}, {CONTROL(sel).filename}); % get selection l = 1;
            % first concatenate the data sets
            tb=[]; fsl = []; fisi = []; ihold=[]; rin=[]; vm = []; nspike=[];
            ks = l*2 - 1; % j - sel(1) + 1;
            l = l + 1;
            s = symbols(ks);
            leg(l,:) = sprintf('%12s',CONTROL(sel2(1)).filename);
            for j = sel(sel2) % just for this one.
                if(~isempty(CONTROL(j).tc_int)) % has analysis been done for this one?
                    if(j == sel(sel2(1)))
                        t0 = min([CONTROL(j).tc_int.time]); % get absolute starting time
                    end;
                    tb = [tb [CONTROL(j).tc_int.time]-t0];
                    fsl = [fsl CONTROL(j).tc_int.fsl];
                    fprintf(1, 'j: %d  #fsl: %d  #tb: %d\n', j, length([CONTROL(j).tc_int.fsl]), ...
                        length([CONTROL(j).tc_int.time]));
                    fisi = [fisi CONTROL(j).tc_int.fisi];
                    ihold = [ihold CONTROL(j).tc_int.ihold];
                    rin = [rin CONTROL(j).tc_int.rin];
                    vm = [vm CONTROL(j).tc_int.vm];
                    nspike = [nspike CONTROL(j).tc_int.nspike];
                end;
            end; % dataset is built

            % FSL FISI plots.

            subplot(hp1);
            plot(tb, fsl, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('FSL (ms)'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);
            legend(leg,2);

            % FISI
            subplot(hp2);
            plot(tb, fisi, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('FISI (ms)'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);


            % Nspikes
            subplot(hp3);
            plot(tb, nspike, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('Nspikes'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);

            % Ihold
            subplot(hp5);
            plot(tb, ihold, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('I_hold (pA)'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);

            % Rin
            subplot(hp6)
            plot(tb, rin, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('Rin (Mohm)'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);
            hold off;

            % Vm
            subplot(hp4);
            plot(tb, vm, 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            hold on;
            ylabel(sprintf('Vm (mV)'), 'Fontsize', 8);
            xlabel('min', 'Fontsize', 8);

        end;


    case 'tc_int_avg'
        % average several tc_ints on one graph.......?
        if(length(sel) < 2)
            return;
        end;
        if(~check_field(CONTROL(sel(1)), 'tc_int') | isempty(CONTROL(sel(1)).tc_int))
            return;
        end;
        l = 1;

        % interpolate the data onto a common time base - use the first set to do it.
        %tc_tb = [CONTROL(sel(1)).tc_int.time];
        tc_tb = [0:(1/12):5-(1/12) 7:(1/6):27-(1/6)]' ; %[CONTROL(sel(1)).FP.time];
        k = 1;
        tfsl = zeros(length(tc_tb), length(sel));
        tfisi = tfsl;
        tvm = tfsl;
        trin = tfsl;
        blx = [];
        for j=sel
            if(check_field(CONTROL(j), 'tc_int') & ~isempty(CONTROL(j).tc_int))
                [tx, bl] = find([CONTROL(j).tc_int.time] < 5 & [CONTROL(j).tc_int.time] >= 0.1);
                fsl = [CONTROL(j).tc_int.fsl];
                fisi = [CONTROL(j).tc_int.fisi];
                vm = [CONTROL(j).tc_int.vm];
                rin = [CONTROL(j).tc_int.rin];
                blfsl = mean(fsl(tx));
                blfisi = mean(fisi(tx));
                blvm = mean(vm(tx));
                blrin = mean(rin(tx));
                tfsl(:,k) = interp1([CONTROL(j).tc_int.time], fsl-blfsl, tc_tb);
                tfisi(:,k) = interp1([CONTROL(j).tc_int.time], fisi-blfisi, tc_tb);
                tvm(:,k) = interp1([CONTROL(j).tc_int.time], vm - blvm, tc_tb);
                trin(:,k) = interp1([CONTROL(j).tc_int.time], rin - blrin, tc_tb);
            end;
            k = k + 1;
        end;
        sf = 1/sqrt(length(sel)); % to get sem.
        [m_fsl, sd_fsl, n_fsl] = mean_nonan(tfsl);
        [m_fisi, sd_fisi, n_fisi] = mean_nonan(tfisi);
        [m_vm, sd_vm, n_vm] = mean_nonan(tvm);
        [m_rin, sd_rin, n_rin] = mean_nonan(trin);

        subplot('Position', [0.1, 0.75, 0.8, 0.2]);
        errorbar(tc_tb, m_fsl, sd_fsl, 'rx')
        ylabel(sprintf('\\Delta FSL\n(mean = %6.2f)', blfsl), 'FontSize', 8);
        subplot('Position', [0.1, 0.5, 0.8, 0.2]);
        errorbar(tc_tb, m_fisi, sd_fisi,'bo')
        ylabel(sprintf('\\Delta FISI\n(mean = %6.2f)', blfisi), 'FontSize', 8);
        subplot('Position', [0.1, 0.3, 0.8, 0.2]);
        errorbar(tc_tb, m_vm, sd_vm,  'ks')
        ylabel(sprintf('\\Delta Vm\n(mean = %6.2f)', blvm), 'FontSize', 8);
        subplot('Position', [0.1, 0.1, 0.8, 0.2]);
        errorbar(tc_tb, m_rin, sd_rin, 'ks')
        ylabel(sprintf('\\Delta Rin\n(mean = %6.2f)', blrin), 'FontSize', 8);

        result = [tc_tb'; m_fsl'; sd_fsl'; n_fsl'; m_fisi'; sd_fisi'; n_fisi'; m_vm' ;sd_vm'; ...
            n_vm'; m_rin'; sd_rin'; n_rin'];
        rtitle = 'tc_tb m_fsl sd_fsl FSL m_fisi sd_fisi FISI m_vm sd_vm Vm m_rin sd_rin Rin';


    case 'FP_TC'
        % multiple FP-TC's on one graph.......?
        l = 1;
        subplot('Position', [0.1, 0.65, 0.8, 0.3]);
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'FP') & ~isempty(CONTROL(j).FP))
                plot([CONTROL(j).FP.time], [CONTROL(j).FP.maxN1], 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
                leg(l,:) = sprintf('%12s',CONTROL(j).filename);
                l = l + 1;
                hold on;
            end;
        end;
        l = 1;
        subplot('Position', [0.1, 0.35, 0.8, 0.3]);
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'FP') & ~isempty(CONTROL(j).FP))
                plot([CONTROL(j).FP.time], [CONTROL(j).FP.maxN2], 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
                leg(l,:) = sprintf('%12s',CONTROL(j).filename);
                l = l + 1;
                hold on;
            end;
        end;
        hold off;
        l = 1;
        subplot('Position', [0.1, 0.05, 0.8, 0.3]);
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'FP') & ~isempty(CONTROL(j).FP))
                plot([CONTROL(j).FP.time], [CONTROL(j).FP.maxP3], 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
                leg(l,:) = sprintf('%12s',CONTROL(j).filename);
                l = l + 1;
                hold on;
            end;
        end;
        hold off;


    case 'FP_TC_avg'
        % average several FP_TCs on one graph.......?
        if(length(sel) < 2)
            return;
        end;
        if(~check_field(CONTROL(sel(1)), 'FP') | isempty(CONTROL(sel(1)).FP))
            return;
        end;
        l = 1;

        % interpolate the data onto a common time base - use the first set to do it.
        tc_tb = [0:(1/3):50-(1/3)] ; %[CONTROL(sel(1)).FP.time];
        tc_tb = [CONTROL(sel(1)).FP.time]/60;
        k = 1;
        tfn2 = zeros(length(sel), length(tc_tb));
        tfp3 = tfn2;
        for j=sel
            if(check_field(CONTROL(j), 'FP') & ~isempty(CONTROL(j).FP))
                tm = [CONTROL(j).FP.time]/60; % convert time to minutes
                if(max(tm) > 50) % if more than 50 minutes .....
                    tzero = 0; % or 20...
                    tb = 0+tzero; te = 30+tzero; to = 20; te2 = 90-tzero;
                else
                    tb = 0; te = 10; to = 0; te2 = 50; % basic test (fp data, from pdac experiments)
                end;
                [tx, bl] = find((tm <= te) & (tm >= tb));
                [tx2] = find((tm <= te2) & (tm >= tb));
                if(min(tx2) > 1)
                    tx2 = [min(tx2)-1; tx2];
                end;
                %         bl
                %         tx
                fn2 = [CONTROL(j).FP.maxN2];
                fp3 = [CONTROL(j).FP.maxP3];
                bln2 = mean(fn2(tx));
                blp3 = mean(fp3(tx));
                %         blfsl = mean([CONTROL(j).tc_int.fsl(bl)]);
                tfn2(k, :) = interp1(tm(tx2)-tb, fn2(tx2)/bln2, tc_tb)';
                tfp3(k, :) = interp1(tm(tx2)-tb, fp3(tx2)/blp3, tc_tb)';
            end;
            k = k + 1;
        end;
        sf = 1/sqrt(length(sel)); % to get sem.
        [m_n2, sd_n2, n_n2] = mean_nonan(tfn2');
        [m_p3, sd_p3, n_p3] = mean_nonan(tfp3');

        subplot('Position', [0.1, 0.65, 0.8, 0.3]);
        he = my_errorbar(tc_tb, m_n2, sd_n2, 'rx', 0.3) ;
        set(he(1), 'MarkerSize', 0.25);
        set(he(2), 'MarkerSize', 3);
        set(get(gca, 'title'), 'string', 'N2');
        set(gca, 'XLim', [0 max(tc_tb)]);
        subplot('Position', [0.1, 0.3, 0.8, 0.3]);
        he=my_errorbar(tc_tb, m_p3, sd_p3, 'bo', 0.3);
        set(he(2), 'MarkerSize', 0.25);
        set(he(1), 'MarkerSize', 3);
        set(gca, 'XLim', [0 max(tc_tb)]);
        set(get(gca, 'title'), 'String', 'P3');
        subplot('Position', [0.1, 0.1, 0.8, 0.2]);
        hp = plot(tc_tb, n_n2, 'rx', tc_tb, n_p3, 'bo');
        set(hp, 'MarkerSize', 1.5);
        set(gca, 'XLim', [0 max(tc_tb)]);
        set(get(gca, 'title'), 'string', 'N2, P3');
        result = [tc_tb'; m_n2'; sd_n2'; n_n2'; m_p3'; sd_p3'; n_p3'];
        rtitle = 'tc_tb m_n2 sd_n2 n_n2 m_p3 sd_pe n_p3';
        return;


    case 'FP_IO' % multiple IO functions, their fits, anad the PPF plots, on one graph.
        hs = cell(2*length(sel), 3);
        jindex = 0;
        ji = 0;
        wavel = [1,3];
        lines = {'-', '--', '.', '-.'};
        colrs = {'k', 'r', 'g', 'b'};
        marks = {'o', 's', '^', '+'};
        for j=sel % copy the information to the main record
            jindex = ji*2+1;
            ji = ji + 1;
            l = 1;
            waven = 0;
            for wave = wavel
                % for the P3, wave = 3, for the n1, wave = 1:
                switch(wave)
                    case 1
                        xpos = 0.05;
                    case 3
                        xpos = 0.40;
                    otherwise
                        xpos = 0.75;
                end;
                if(isempty(hs{jindex, wave}))
                    hs{jindex, wave} = subplot('Position', [xpos, 0.52, 0.3, 0.3]);
                    hold on;
                else
                    subplot(hs{jindex, wave});
                    hold on;
                end;
                ks = l; % (j - sel(1)) - 1;
                s = symbols(ks);
                if(check_field(CONTROL(j), 'FP') & ~isempty(CONTROL(j).FP))
                    x=CONTROL(j).FP.S1;
                    %  x=x(1);
                    hp=semilogx(CONTROL(j).FP.IO, x(:,wave));
                    hold on;
                    set(hp, 'Color', colrs{ji}, ...
                        'Marker', marks{wave}, 'MarkerFaceColor', colrs{ji}, 'markersize', 3, 'LineStyle', 'none');
                    hp=semilogx(CONTROL(j).FP.S1_fitx(wave,:),CONTROL(j).FP.S1_fity2(wave,:));
                    set(hp, 'Color', colrs{ji}, ...
                        'Marker', 'none', 'MarkerFaceColor', colrs{ji}, 'LineStyle', '-')
                    leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                    sn = symbols(ks+1);
                    y=CONTROL(j).FP.S2;
                    %  y=y(1);
                    hp=semilogx(CONTROL(j).FP.IO, y(:,wave));
                    set(hp, 'Color', colrs{ji}, ...
                        'Marker', marks{wave}, 'MarkerFaceColor', 'none', 'markersize', 3, 'LineStyle', 'none');
                    hp=semilogx(CONTROL(j).FP.S2_fitx(wave,:),CONTROL(j).FP.S2_fity2(wave,:));
                    set(hp, 'Color', colrs{ji}, ...
                        'Marker', 'none', 'MarkerFaceColor', colrs{ji}, 'LineStyle', '-')
                    leg(l+1,:) = sprintf('%12s',CONTROL(j).reclist);
                    set(gca, 'Xscale', 'log');
                    set(gca, 'fontsize', 8);
                    ht = get(gca, 'title');
                    switch(wave)
                        case 1
                            set(ht, 'string', 'N1');
                        case 2
                            set(ht, 'string', 'N2');
                        case 3
                            set(ht, 'string', 'P3');
                        otherwise
                    end;

                end;
                % legend(leg,1);
                if(isempty(hs{jindex+1, wave}))
                    hs{jindex+1, wave} = subplot('Position', [xpos, 0.05, 0.3, 0.3]);
                    hold on;
                else
                    subplot(hs{jindex+1, wave});
                    hold on;
                end;

                %         subplot('Position', [xpos, 0.05, 0.6, 0.44]);
                PP = zeros(size(x,1),1);
                for k = 1:size(x, 1)
                    if(y(k,3) < 2*CONTROL(j).FP.stdn | x(k,wave) < 2*CONTROL(j).FP.stdn)
                        PP(k) = NaN;
                    else
                        PP(k) = y(k,wave)./x(k,wave);
                    end;
                end;
                hp=semilogx(CONTROL(j).FP.IO, PP);
                hold on;
                set(hp, 'Color', colrs{ji}, ...
                    'Marker', marks{wave}, 'MarkerFaceColor', colrs{ji}, 'markersize', 2, 'LineStyle', 'none');
                hp=semilogx(CONTROL(j).FP.S1_fitx(wave,:),CONTROL(j).FP.S2_fity2(wave,:)./CONTROL(j).FP.S1_fity2(wave,:));
                set(hp, ...
                    'Color', colrs{ji}, ...
                    'Marker', 'none', 'MarkerFaceColor', colrs{ji}, 'LineStyle', '-');
                leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                u=get(gca, 'Ylim');
                set(gca, 'Ylim', [-0.5 u(2)]);
                l = l + 1;
                set(gca, 'Xscale', 'log');
                set(gca, 'fontsize', 8);
                hold off;
                ht = get(gca, 'title');
                switch(wave)
                    case 1
                        set(ht, 'string', 'N1_p_p');
                    case 2
                        set(ht, 'string', 'N2_p_p');
                    case 3
                        set(ht, 'string', 'P3_p_p');
                    otherwise
                end;
            end;


            %legend(leg,2);

        end;
        % now take the second selection and ratio it against the first, as a function of IO
        if(length(sel) == 2) % only if length is 2.
            subplot('Position', [0.75, 0.52, 0.25, 0.35]);
            j1 = sel(1);
            j2 = sel(2);
            x=CONTROL(j1).FP.S1;
            y = CONTROL(j2).FP.S1;
            semilogx(CONTROL(j1).FP.IO, y(:, 3)./x(:,3), 'color', 'k', 'linestyle', 'none', 'marker', '^', 'markersize', 2);
            hold on;
            hp=semilogx(CONTROL(j1).FP.S1_fitx(3,:), CONTROL(j2).FP.S1_fity2(3,:)./CONTROL(j1).FP.S1_fity2(3,:));
            set(hp, 'color', 'r', 'linestyle', '-', 'marker', 'none');
            set(gca, 'Ylim', [0 u(2)]);
            set(gca, 'fontsize', 8);
            hold off;
            ht = get(gca, 'title');
            set(ht, 'string', 'P3_{pp} - Ctl vs PDAc');

            %	subplot('Position', [0.70, 0.05, 0.28, 0.4]);
            %   j1 = sel(1);
            %   j2 = sel(2);
            %  x = CONTROL(j1).FP.S1;
            %  y = CONTROL(j2).FP.S2;
            %   plot(x(:, 3), y(:, 3)./x(:,3), 'color', 'k', 'linestyle', 'none', 'marker', 'o', 'markersize', 2); % ppf vs r1
            %   hold on;
            %   plot(y(:, 3), y(:, 3)./x(:,3), 'color', 'r', 'linestyle', 'none', 'marker', '^', ...
            %      'markerfacecolor', 'none', 'markersize', 2); % ppf vs r2
            %  hold off;

        end;

        subplot('Position',[0.07,0.95,0.83,0.05])
        axis([0,1,0,1])
        axis('off')
        text(0.5,0.1,sprintf('%-12s R[%s,%s] %-8s', ...
            CONTROL(sel(1)).filename, CONTROL(sel(1)).reclist, CONTROL(sel(2)).reclist, CONTROL(sel(1)).protocol), ...
            'Fontsize', 8, ...
            'horizontalalignment', 'center');


        % allow us to store the results.
        %result = [tc_tb'; m_n2; sd_n2; n_n2; m_p3; sd_p3; n_p3];
        %rtitle = 'tc_tb m_n2 sd_n2 n_n2 m_p3 sd_pe n_p3';

        % EPSC vs time in a train from EPSC_depression protocol
    case('ESPC_Dep') % multiply selected measurements in one graph (typo in mainscreen for funny name)
        l = 1;
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            if(check_field(CONTROL(j), 'EPSC_Dep') & ~isempty(CONTROL(j).EPSC_Dep))
                n=[1:CONTROL(j).EPSC_Dep.npsc];
                xp = CONTROL(j).stim_time;
                yp = -CONTROL(j).EPSC_Dep.Im;
                ye = sqrt(CONTROL(j).EPSC_Dep.Ivar);
                ipi = mean(diff(CONTROL(j).stim_time));
                first = 3;
                xfit = [xp(first):xp(end)];
                hx = errorbar(xp, yp, ye);
                set(hx, 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'LineStyle', '-');
                leg{l} = sprintf('%12s',CONTROL(j).reclist);
                %  [amp0 amp1 tau] = expfit(xp(first:end), yp(first:end));
                amp0 = 1000;
                amp1 = 20000;
                tau = 25;
                lb = [0 0 0];
                ub = [amp0*20 amp1*50 tau*20];
                if(amp0 < 0)
                    amp0 = 0;
                    ub(1) = 20000;
                end;
                if(amp1 < 0)
                    amp1 = 1000;
                    ub(2) = 20000;
                end;
                % use that to seed the m-l fit
                [fpar, chisq, niter, epsc_fit1] = mrqfit('exponential', [amp0 amp1 tau], xp(first:end), yp(first:end), ye(first:end), ...
                    [], lb, ub, 500, []);
                if(isempty(epsc_fit1))
                    epsc_fit = amp0 + amp1*exp(-xfit(2:end)/tau);
                end
                fp1=fpar;
                cerr1=chisq;
                amp10 = fpar(1); amp11=fpar(2); tau11 = fpar(3);
                ub = [amp10*20 amp11*200 tau11*10 amp11*5000 tau11*10];
                lb = [0 0 0 0 0]; % require all to be positive
                [fpar2, chisq2, niter2, epsc_fit2] = mrqfit('exponential', [amp10 amp11/2 tau11*2 amp11/2 tau11/4], xp(first:end), yp(first:end), ye(first:end), ...
                    [], lb, ub, 500, []);
                amp2 = fpar2(1); amp21=fpar2(2); tau21 = fpar2(3); amp22 = fpar2(4); tau22 = fpar2(5);
                %time_fit = time_fit+tdel+deadwin;
                if(tau21 > 50 | tau22 > 50)
                    chisq2=5*chisq; % make it a worse fit.
                end
                fp2=fpar2;
                cerr2=chisq2;
                l = l + 1;
                hold on;
                %  if(mean(cerr2./cerr1) < 1 & amp22 > 0)  % require better fit and positive amplitude
                yfit2 = amp2 + amp21*exp(-xfit/tau21) + amp22*exp(-xfit/tau22);
                plot(xfit, yfit2, 'r-');
                %  else
                yfit1 = amp10 + amp11*exp(-xfit/tau11);
                plot(xfit, yfit1, 'c-'); % double exponential
                yfit0 = amp0 + amp1*exp(-xfit/tau);
                % plot(xfit, yfit0, 'k-'); % initial cheb fit
                %  end;
                %yfit = amp0 + amp1*exp(-xfit/tau);
            else
                QueMessage('Empty EPSC_dep in file?');
            end;
        end;
        hold off;
        if(~isempty(leg))
            legend(leg,2, 'Location', 'northeast');
        end;

    case 'EPSP_TC_OLD'
        for j = sel
            if(~check_field(CONTROL(j), 'PSP_TC') | isempty(CONTROL(j).PSP_TC))
                QueMessage('Compare_plot: no EPSP_TC data in field?');
                return;
            end;
        end;

        % EPSP_TC's - are combined in time...
        l = 1;
        wrapflag = 0;
        subplot('Position', [0.1, 0.72, 0.8, 0.25]);
        j = sel(1);
        zt = CONTROL(j).PSP_TC.zt;
        ts = zt(1)/60; % convert to minutes...
        t0 = 0;
        ib = find([CONTROL(j).PSP_TC.T_epsp_1]+t0 < 10 & [CONTROL(j).PSP_TC.T_epsp_1]+t0 >= 5);
        vepsp = [];
        tc = [];
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            zt = CONTROL(j).PSP_TC.zt;
            if(wrapflag)
                zt = zt + 60 * 60 * 24;
            end;
            x=find(diff(zt)<0);
            if(~isempty(x)) % wrap around......
                zt(x:end) = zt(x:end)+60*60*24; % add a day
                fprintf('wrapflag set\n');
                wrapflag = 1;
            end;
            if(j > sel(1))
                t0 = zt(1)/60 - ts; % convert to minutes...
            end;
            plot([CONTROL(j).PSP_TC.T_epsp_1]+t0, [CONTROL(j).PSP_TC.V_epsp_1], 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            vepsp = [vepsp CONTROL(j).PSP_TC.V_epsp_1];
            tc = [tc CONTROL(j).PSP_TC.T_epsp_1+t0];
            leg(l,:) = sprintf('%12s',CONTROL(j).filename);
            l = l + 1;
            hold on;
        end;
        [vme, vvar] = mean_var([CONTROL(sel(1)).PSP_TC.V_epsp_1(ib)]);
        vvar = sqrt(vvar);
        u = get(gca, 'XLim');
        plot(u, [vme vme], 'linestyle', '-', 'color', 'k');
        plot(u, [vme+vvar vme+vvar], 'linestyle', '--', 'color', 'k');
        plot(u, [vme-vvar vme-vvar], 'linestyle', '--', 'color', 'k');
        u = get(gca, 'Ylim');
        set(gca, 'Ylim', [0 u(2)]);
        hold off;
        ylabel(sprintf('V_{EPSP} (mV)'), 'Fontsize', 8);
        xlabel('min', 'Fontsize', 8);
        legend(leg,2);
        % generate normalized data for "result"
        vepsp = 100*vepsp/vme; % express as percent control

        subplot('Position', [0.1, 0.45, 0.8, 0.25]);
        t0 = 0;
        wrapflag = 0;
        vsepsp = [];
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            zt = CONTROL(j).PSP_TC.zt;
            if(wrapflag)
                zt = zt + 60 * 60 * 24;
            end;
            x=find(diff(zt)<0);
            if(~isempty(x)) % wrap around......
                zt(x:end) = zt(x:end)+60*60*24; % add a day
                fprintf('wrapflag set\n');
                wrapflag = 1;
            end;
            if(j > sel(1))
                t0 = zt(1)/60 - ts; % convert to minutes...
            end;
            plot([CONTROL(j).PSP_TC.T_s_epsp_1]+t0, [CONTROL(j).PSP_TC.S_epsp_1], 'Color', s.c, ...
                'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
            %         leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
            vsepsp = [vsepsp CONTROL(j).PSP_TC.S_epsp_1];
            l = l + 1;
            hold on;
        end;
        [vms, vvar] = mean_var([CONTROL(sel(1)).PSP_TC.S_epsp_1(ib)]);
        vvar = sqrt(vvar);
        u = get(gca, 'XLim');
        plot(u, [vms vms], 'linestyle', '-', 'color', 'k');
        plot(u, [vms+vvar vms+vvar], 'linestyle', '--', 'color', 'k');
        plot(u, [vms-vvar vms-vvar], 'linestyle', '--', 'color', 'k');
        hold off;
        vsepsp = vsepsp*100/vms;

        u = get(gca, 'Ylim');
        set(gca, 'Ylim', [0 u(2)]);
        ylabel(sprintf('EPSP Slope (mV/ms)'), 'Fontsize', 8);
        xlabel('min', 'Fontsize', 8);
        hold off;

        vi = [];
        bl = 1;
        t0 = 0;
        wrapflag = 0;
        subplot('Position', [0.1, 0.18, 0.8, 0.25]);
        vipsp = [];
        for j=sel % copy the information to the main record
            ks = l*2 - 1; % j - sel(1) + 1;
            s = symbols(ks);
            zt = CONTROL(j).PSP_TC.zt;
            if(wrapflag)
                zt = zt + 60 * 60 * 24;
            end;
            x=find(diff(zt)<0);
            if(~isempty(x)) % wrap around......
                zt(x:end) = zt(x:end)+60*60*24; % add a day
                fprintf('wrapflag set\n');
                wrapflag = 1;
            end;
            if(j > sel(1))
                t0 = zt(1)/60 - ts; % convert to minutes...
            end;
            if(check_field(CONTROL(j), 'PSP_TC') & ~isempty(CONTROL(j).PSP_TC) & length(CONTROL(j).PSP_TC.V_ipsp_1) > 1)
                % special for IPSPs: only plot if an IPSP is detected; e.g., if V_ipsp_1 < 0.
                ipsp = [CONTROL(j).PSP_TC.V_ipsp_1];
                ipsp(find(ipsp < 0)) = NaN;
                plot([CONTROL(j).PSP_TC.T_s_epsp_1]+t0, ipsp, 'Color', s.c, ...
                    'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
                vipsp = -[vipsp CONTROL(j).PSP_TC.V_ipsp_1];
                %      plot(ta, va, 'Color', s.c, ...
                %         'Marker', s.s, 'MarkerFaceColor', s.f, 'MarkerSize', 1.5, 'LineStyle', '-');
                %        leg(l,:) = sprintf('%12s',CONTROL(j).reclist);
                l = l + 1;
                hold on;
            end;
        end;
        % ib = find(ta0 < 10 & ta0 >= 5);
        ipm = [CONTROL(sel(1)).PSP_TC.V_ipsp_1(ib)];
        ipm(find(ipm > 0)) = NaN;
        if(any(ipm ~= NaN))
            [vmi, vvar] = mean_var(ipm);
            vvar = sqrt(vvar);
            u = get(gca, 'XLim');
            plot(u, -[vmi vmi], 'linestyle', '-', 'color', 'k');
            plot(u, -[vmi+vvar vmi+vvar], 'linestyle', '--', 'color', 'k');
            plot(u, -[vmi-vvar vmi-vvar], 'linestyle', '--', 'color', 'k');
            u = get(gca, 'Ylim');
            if(u(1) == u(2))
                u(2) = u(1)+1;
            end;
            set(gca, 'Ylim', [0 u(2)]);
        end;
        hold off;
        vipsp = 100*vipsp/vmi;

        ylabel(sprintf('V_{IPSP} (mV)'), 'Fontsize', 8);
        xlabel('min', 'Fontsize', 8);
        hold off;
        result = [tc;vepsp;vsepsp;vipsp];
        rtitle = 'Time VEPSP VSlope VIPSP';

        % new version with block averaging...
    case 'EPSP_TC'
        for j = sel
            if(~check_field(CONTROL(j), 'PSP_TC') | isempty(CONTROL(j).PSP_TC))
                QueMessage('Compare_plot: no EPSP_TC data in field?');
                return;
            end;
        end;
        if(arg2 > 0) % skip the dialog box
            % bring up dialog box.
                clear_flag = 0;
                accum_flag = 0;
                block_window = 1;
                tbeg = 0;
                tend = 5;
                [result, rtitle] = EPSC_TCP(clear_flag, accum_flag, block_window, tbeg, tend, result, rtitle, state_arg);
                return;
        end;
        
        hctc= findobj('tag', 'combine_epsptc');
        if(isempty(hctc) | ~ishandle(hctc))
            open('combine.fig');
            hctc= findobj('tag', 'combine_epsptc');
        end;
        if(nargin > 1)
            state_arg = varargin{2}; % state_arg holds calling state - 0 means use voltage out, 1 means use slope
        else
            state_arg = 0;
        end;
        EPSC_STATE=state_arg; % transfer to global variable

        figure(hctc); % visible.
        return; % wait for command...

    case 'EPSP_TC_GO' % ok to proceed.
        % Multiple EPSC_TC's are combined as a time sequence
        % and plotted. Each time point is reduced to an average.
        % from surrounding points.
        %
        % read parameters
        block_window = get_window_arg('combine_window', 1); % block window in time units (minutes)
        tbeg = get_window_arg('combine_begin', 0);
        tend = get_window_arg('combine_end', 5);
        accum_flag = get_window_arg('combine_acc', 0, 1); % read accumulate checkbox
        clear_flag = get_window_arg('combine_clear', 0, 1); % read clear flag
        hctc= findobj('tag', 'combine_epsptc');
        close(hctc);
        [result, rtitle] = EPSC_TCP(clear_flag, accum_flag, block_window, tbeg, tend, result, rtitle, state_arg);

    case 'EPSC_TC_AUTO'
        if(nargin > 1)
            state_arg = varargin{2}; % state_arg holds calling state - 0 means use voltage out, 1 means use slope
        else
            state_arg = 0;
        end;
        EPSC_STATE=state_arg; % transfer to global variable
        % [result, rtitle] = EPSC_TCP(0, 1, 0.16666666, 0, 10, result, rtitle, state_arg);
        [result, rtitle] = EPSC_TCP(0, 1, 1, 3, 10, result, rtitle, state_arg);
        TRES = result;

    case 'export'
        %result
        if(isempty(result))
            msgbox('No result?');
            return;
        end;
        switch computer
            case {'MAC', 'MACI'}
                opath = '/users/pbmanis/Desktop/';
                % PC:
            case 'PCWIN'
                opath = 'c:\data\extracted\';
            otherwise % assume glnx86 ...
                opath = '/home/pmanis/extracteddata/';
        end;
        rlist = seq_parse(CONTROL(sel(1)).reclist)
        if(length(rtitle) < 12)
            fn = sprintf('%s%s_%d_%s.dat', opath, CONTROL(sel(1)).filename,  rlist{1}(1), rtitle)
        else
            fn = sprintf('%s%s_%d.dat', opath, CONTROL(sel(1)).filename,  rlist{1}(1))
        end;
        hw=fopen(fn, 'w');
        if(isempty(hw) | hw == [])
            fprintf(1, 'Unable to open output file');
            return;
        end;
        rl = size(result);
        fmtstring = [];
        if(~iscell(result))
            for i = 1:rl(1)
                fmtstring = [fmtstring '%f, '];
            end;
        else
            for i = 1:length(result)
                fmtstring = [fmtstring '%f, %f, '];
            end;
        end;
        fmtstring = [fmtstring '\r\n'];
        fid = fopen(fn, 'w');
        fprintf(fid, '%s \r\n', rtitle);
        if(iscell(result))
            done = 0;
            row = 1;
            while(~done)
                pflg = 0;
                for i = 1:length(result)
                    a=[result{i}]; % get the data itself
                    if(size(a,2) >= row)
                        fprintf(fid, '%f, %f, ', a(:,row)');
                        pflg = 1; % signal we printed at least one thing...
                    else
                        fprintf(fid, ', , '); % just the placeholders, no data
                    end;
                end;
                fprintf(fid, '\r\n'); % end the row
                row = row + 1; % move to next row
                if(pflg== 0) done = 1; end;
            end;

        else
            for i = 1:rl(2)
                fprintf(fid, fmtstring, result(:,i));
            end;
        end;
        fclose(fid);
    otherwise
end

return;

function yp = nextln(yp)
yp = yp - 0.05;
if(yp < 0)
    yp = 0;
elseif(yp > 1)
    yp = 1;
end;


    function [tout] = window_average(t, wsize)
        %
        % compute block averages of v data over wsize window (seconds)
        % using window size wsize. The data are returned in new arrays
        % t
        %
        nd=length(t);
        j=1;
        for i = 1:wsize:nd
            tout(j) = mean(t(i:i+wsize-1));
            j = j+1;
        end;


        function [hax] = plot_witherror(t, v, verr, ib, title,xtick)
            % plot data t and v, with standard error bars verr
            % also insert baseline error measurement from the points defined by indices ib
            % as line and +/- limits.
            if(any(v ~= NaN))
                [vm, vvar] = mean_var(v(ib));
                vvar = sqrt(vvar);

                if(any(verr) > 0)
                    hplot = my_errorbar(t, v, verr, verr, 'ko', 0.025);
                end;
                hold on;
                hplot = plot(t, v, 'ko');
                set(hplot, 'Markersize', 1.5, 'MarkerFaceColor', 'k', 'markeredgecolor', 'k');
                hold on;
                %set(hplot, 'color', 'black');
                %plot(t, v, 'Color', 'black', ...
                %'Marker', 'o', 'MarkerFaceColor', 'k', 'MarkerSize', 1.5, 'LineStyle', 'none')
                set(gca, 'Xlim', [-1 50]);
                u = get(gca, 'XLim');
                plot(u, [vm vm], 'linestyle', '-', 'color', 'k');
                plot(u, [vm+vvar vm+vvar], 'linestyle', '--', 'color', 'k');
                plot(u, [vm-vvar vm-vvar], 'linestyle', '--', 'color', 'k');
                % u = get(gca, 'Ylim');
                % set(gca, 'Xlim', [-1 50]);
                % set(gca, 'Ylim', [u(1) u(2)]);
                hold off;
                ylabel(sprintf(title), 'Fontsize', 8, 'interpreter', 'tex');
                xlabel('time (min)', 'Fontsize', 8);
                hax = gca; % return handle to axes
                if(nargin < 5)
                    set(gca,'XTick',[ ]);
                end;
                return
            end;


            function [result, rtitle] = EPSC_TCP(clear_flag, accum_flag, block_window, tbeg, tend, result, rtitle, state_arg)

                global CONTROL

                msize = 3;
                fsize = 7;

                h2=findobj('Tag', 'CTL'); % first get the selection number
                sel = get(h2, 'Value');

                if(clear_flag | ~accum_flag)
                    result = [];
                end;
                l = 1;
                wrapflag = 0;
                j = sel(1);
                zt = CONTROL(j).PSP_TC.zt;
                ts = zt(1)/60; % convert to minutes...
                t0 = 0;
                tc = [];
                tc2 = [];
                % analyzed the EPSP/IPSP
                vepsp = [];
                vsepsp = [];
                vipsp = [];
                % and the baseline information about stability of the cell
                rin = [];
                tau = [];
                vm = [];
                ih = [];
                t0 = CONTROL(sel(1)).PSP_TC.zt; % note: ftime is the ENDING time, not start time for the data
                dmode = CONTROL(sel(1)).mode; % get the acquisition mode
                t0 = t0(1);
                for j=sel %
                    if(strmatch(CONTROL(j).protocol, {'timing_base', 'hpc_ltpvc_base'}))
                        ks = l*2 - 1; % j - sel(1) + 1;
                        s = symbols(ks);
                        %         zt = CONTROL(j).PSP_TC.zt;
                        %         if(wrapflag)
                        %             zt = zt + 60 * 60 * 24;
                        %         end;
                        %         x=find(diff(zt)<0);
                        %         if(~isempty(x)) % wrap around......
                        %             zt(x:end) = zt(x:end)+60*60*24; % add a day
                        %             fprintf('wrapflag set\n');
                        %             wrapflag = 1;
                        %         end;
                        %         if(j > sel(1))
                        %             t0 = zt(1)/60 - ts; % convert to minutes...
                        %         end;
                        vepsp = [vepsp CONTROL(j).PSP_TC.V_epsp_1'];
                        vsepsp = [vsepsp CONTROL(j).PSP_TC.S_epsp_1'];
                        vipsp = -[vipsp CONTROL(j).PSP_TC.V_ipsp_1'];
                        if(check_field(CONTROL(j).PSP_TC, 'Rin'))
                            rin = [rin CONTROL(j).PSP_TC.Rin];
                            tau = [tau CONTROL(j).PSP_TC.Tau];
                            vm = [vm CONTROL(j).PSP_TC.Vm];
                            ih = [ih CONTROL(j).PSP_TC.Ih];
                        else
                            dlen = length(CONTROL(j).PSP_TC.V_epsp_1);
                            padd = zeros(1, dlen);
                            rin = [rin padd];
                            tau = [tau padd];
                            vm = [vm padd];
                            ih = [ih padd];
                        end;
                        %  tc = [tc CONTROL(j).PSP_TC.T_epsp_1+t0];
                        %tc2 = [tc2 CONTROL(j).PSP_TC.T+(CONTROL(j).PSP_TC.ftime-max(CONTROL(j).PSP_TC.T)-t0)]; % time base for these measurements
                        tc2 = [tc2 CONTROL(j).PSP_TC.zt-t0];
                        l = l + 1;
                        hold on;
                    end;
                end;
                tc2 = tc2/60;

                sp1=subplot('Position', [0.1, 0.94, 0.8, 0.05]);
                axis([0,1,0,1])
                axis('off')
                sf = sel(1);
                txtstring = sprintf('%-12s    R: %d    %-20s',CONTROL(sf).filename, ...
                    CONTROL(sf).recbeg, CONTROL(sf).protocol);
                ht=findobj('tag', 'compare_text');
                if(~isempty(ht))
                    delete(ht);
                end;

                ht=text(0,0.5, txtstring, 'Fontsize', 8, 'interpreter', 'none', 'tag', 'compare_text');

                subplot('Position', [0.1, 0.9, 0.69, 0.039]);
                axis([0,1,0,1])
                axis('off')
                sf = sel(1);
                if(~check_field(CONTROL(sf),'File_Comment') | isempty(CONTROL(sf).File_Comment))
                    truefile = [CONTROL(sf).path CONTROL(sf).filename];
                    [DFILE, data, err] = datac2mat(truefile, []); % just read the header - empty RL list
                    CONTROL(sf).File_Comment = DFILE.comment;
                end;

                com = wrapit(CONTROL(sf).File_Comment); % create a wrapped string.
                txtcom = sprintf('%s', com);
                ht=findobj('tag', 'compare_com');
                if(~isempty(ht))
                    delete(ht);
                end;
                ht=text(0, 0.5, txtcom,  'Fontsize', 8, 'interpreter', 'none', 'tag', 'compare_com');

                [tc_b, vep_b, vep_std, vep_n] = block_average(tc2, vepsp, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);
                subplot('Position', [0.1, 0.65, 0.8, 0.25]);
                if(strcmpi(dmode, 'cc'))
                    lab1 = 'V_{EPSP} (mV)';
                    lab2 = 'EPSP Slope (mV/ms)';
                else
                    lab1 = 'I_{EPSC} (pA)';
                    lab2 = 'EPSC slope (pA/ms)';
                end;
                plot_witherror(tc_b, vep_b, vep_std, ib, lab1);
                u = get(gca, 'Ylim');
                set(gca, 'Ylim', [0 u(2)]);
                set(gca, 'fontsize', fsize);
                set(gca, 'box', 'off', 'tickdir', 'out', 'Xlabel', text('string', ' '));
                % generate normalized data for "result"
                vep_b = 100*vep_b/mean_var(vep_b(ib)); % express as percent control

                [tc_b, vsep_b, vsep_std, vsep_n] = block_average(tc2, vsepsp, block_window);
                subplot('Position', [0.1, 0.35, 0.8, 0.25]);
                plot_witherror(tc_b, vsep_b, vsep_std, ib, lab2);
                u = get(gca, 'Ylim');
                set(gca, 'Ylim', [0 u(2)]);
                set(gca, 'fontsize', fsize);
                set(gca, 'box', 'off', 'tickdir', 'out', 'Xticklabel', '', 'Xlabel', text('string', ' '));
                % generate normalized data for "result"
                vsep_b = 100*vsep_b/mean_var(vsep_b(ib)); % express as percent control

                [tc_b, vip_b, vip_std, vip_n] = block_average(tc2, vipsp, block_window);
                ib = find(tc_b < 10 & tc_b >= 5);
                %subplot('Position', [0.1, 0.38, 0.8, 0.2]);
                vip_b(find(vip_b < 0)) = NaN;
                vip_std(find(vip_b < 0)) = NaN;
                %hax = plot_witherror(tc_b, vip_b, vip_std, ib, 'V_{IPSP} (mV)');


                [tc_b, vm_b, vm_std, vm_n] = block_average(tc2, vm, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);
                [tc_b, ih_b, ih_std, ih_n] = block_average(tc2, ih, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);

                subplot('Position', [0.1, 0.24, 0.8, 0.08]);
                if(strcmpi('cc', dmode))

                    plot_witherror(tc_b, vm_b, vm_std, ib, 'V_{m} (mV)');
                else
                    plot_witherror(tc_b, ih_b, ih_std, ib, 'I_{hold} (pA)');
                end;
                % generate normalized data for "result"
                u=get(gca, 'YLim');
                dy = abs(diff(u));
                if(dy < 15)
                    um=(u(1)+u(2))/2;
                    set(gca, 'YLim', [um-7.5 um+7.5]);
                end;
                set(gca, 'fontsize', fsize);
                set(gca, 'box', 'off', 'tickdir', 'out', 'XTickLabel', '');
                vm_b = 100*vm_b/mean_var(vm_b(ib)); % express as percent control

                [tc_b, ih_b, ih_std, ih_n] = block_average(tc2, ih, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);
                %subplot('Position', [0.1, 0.21, 0.8, 0.07]);
                %plot_witherror(tc_b, ih_b, ih_std, ib, 'I_{H} (pA)');
                %set(gca, 'fontsize', fsize);
                %set(gca, 'box', 'off', 'tickdir', 'out', 'XTickLabel', '')
                % generate normalized data for "result"
                ih_b = 100*ih_b/mean_var(ih_b(ib)); % express as percent control

                [tc_b, rin_b, rin_std, rin_n] = block_average(tc2, rin, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);

                subplot('Position', [0.1, 0.14, 0.8, 0.08]);
                plot_witherror(tc_b, rin_b, rin_std, ib, 'R_{in} (Mohm)');
                u=get(gca, 'YLim');
                dy = abs(diff(u));
                if(dy < 50)
                    um=(u(1)+u(2))/2;
                    set(gca, 'YLim', [um-25 um+25]);
                end;
                set(gca, 'fontsize', fsize);
                set(gca, 'box', 'off', 'tickdir', 'out')
                % generate normalized data for "result"
                rin_b = 100*rin_b/mean_var(rin_b(ib)); % express as percent control

                [tc_b, tau_b, tau_std, tau_n] = block_average(tc2, tau, block_window);
                ib = find(tc_b < tend & tc_b >= tbeg);
                %subplot('Position', [0.1, 0.07, 0.8, 0.07]);
                %plot_witherror(tc_b, tau_b, tau_std, ib, 'Tau_m (msec)', 1);
                %set(gca, 'fontsize', fsize);
                %set(gca, 'box', 'off', 'tickdir', 'out')
                % generate normalized data for "result"
                tau_b = 100*tau_b/mean_var(tau_b(ib)); % express as percent control

                % generate normalized data for "result"
                if(any(vip_b(ib) ~= NaN) & all(vip_b(ib) > 0))
                    vip_b = 100*vip_b/mean_var(vip_b(ib)); % express as percent control
                else
                    vip_b = zeros(1, length(vip_b));
                end;

                if(~isempty(result))
                    if(length(result) > length(vep_b))
                        vep_b(length(vep_b)+1:length(result)) = NaN;
                        vsep_b(length(vsep_b)+1:length(result)) = NaN;
                        tc_b(length(tc_b)+1:length(result)) = NaN;
                    else
                        vep_b = vep_b(1:length(result));
                        vsep_b = vsep_b(1:length(result));
                        tc_b = tc_b(1:length(result));
                    end;
                end;
                % fprintf(1, '%d  l(r) = %d\n', length(vep_b), length(result))
                if(length(tc_b) > length(vsep_b))
                    tc_b = tc_b(1:length(vsep_b));
                end;

                if(isempty(result) | clear_flag)
                    switch(state_arg)
                        case 0
                            result = [tc_b;vep_b];
                        case 1
                            result = [tc_b; vsep_b];
                        otherwise
                            result = [tc_b;vep_b];
                    end;

                    rtitle = [sprintf('Time %s ', CONTROL(sel(1)).filename)];
                else
                    switch(state_arg)
                        case 0
                            result = [result; vep_b];
                        case 1
                            result = [result; vsep_b];
                        otherwise
                            result = [result; vep_b];
                    end;
                    rtitle = [rtitle sprintf('%s ', CONTROL(sel(1)).filename)];
                end;
                orient landscape
                pflag = getplotflag;
                if (pflag > 0)
                    print -dljet3;
                end

                function [y] = gaussfunc(x, fp)
                    %
                    % calculate a gaussian based on x and FP
                    %
                    y = fp(1) + (fp(2)/(fp(4)*sqrt(2*pi)))*exp(-((x-fp(3)).^2)/(2*fp(4)^2));


                    function stro = wrapit(stri)
                        % very simple text wrapping function - splits words....
                        k = 1;
                        for i = 1:length(stri)
                            if(mod(i, 80) ~= 0)
                                stro(k) = stri(i);
                            else
                                stro(k) = [10];
                                k = k + 1;
                                stro(k) = stri(i);
                            end;
                            k = k + 1;
                        end;
